#pragma once

#include "Matrix.h"
#include <vector>
#include <deque>

#define COMPLETE_SPLINE 0 // 完全样条
#define NOTAKNOT_SPLINE 1 // 无结样条
#define PERIODIC_SPLINE 2 // 周期样条

/**
 * @file Polynomial.h
 * @author xingyifan
 * @date 2022-12-09 10:37
 *
 * @description: 一元多项式
 */

namespace UniPolynomial
{
    using BCType = int;
    using Table = std::vector<std::vector<Real>>;

    /**
     * @file Polynomial.h
     * @author xingyifan
     * @date 2022-11-30 16:04
     *
     * @description: 多项式类型
     */

    class Polynomial
    {
    protected:
        std::vector<Real> m_coefficients; // 系数向量

    public:
        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-07-26 13:34
         *
         * @description: 构造器，允许阶数为 0
         */

        explicit Polynomial(int o = 0);            // 限定显式调用，防止 operator* 中的隐式转换
        Polynomial(std::initializer_list<Real> l); // 初始化列表构造
        Polynomial(std::vector<Real> &v);          // 向量构造
        Polynomial(const Polynomial &poly);        // 拷贝构造

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-07-26 13:35
         *
         * @description: 重载操作符
         */

        // 多项式赋值
        void operator=(const Polynomial &poly);

        // 多项式加减法
        Polynomial operator+(const Polynomial &poly) const;
        Polynomial operator-(const Polynomial &poly) const;

        // 多项式乘常数
        Polynomial operator*(Real k) const;
        Polynomial operator/(Real k) const;

        Polynomial operator*(const Polynomial &poly) const; // 多项式相乘
        Polynomial operator/(const Polynomial &poly) const; // 带余除法，获得商多项式

        Real operator()(Real x) const; // 多项式求值

        // 输出函数，此函数就放在这里
        friend std::ostream &operator<<(std::ostream &o, const Polynomial &poly)
        {
            for (int i = 0; i < poly.order(); i++)
                o << poly.coeff(i) << "x^" << i << " ";
            return o;
        }

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-07-27 15:46
         *
         * @description: 多项式操作
         */

        void push(Real k);                    // 推入高阶项
        void reduce();                        // 多项式降阶，将系数为 0 的高阶项去除
        Real intergral(Real a, Real b) const; // 定积分
        Polynomial intergral() const;         // 不定积分
        Polynomial derivative() const;        // 导函数

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-07-26 13:36
         *
         * @description: 获得多项式信息
         */

        int order() const;       // 获取多项式的阶
        Real coeff(int n) const; // 获得对应次数的系数

        // 多项式的拉普拉斯阵，为了方便起见，这里将矩阵取负，这样通过幂法解出的最大特征值就是多项式的根
        LinearAlgebra::Matrix Laplace() const;
    };

    /**
     * @file Polynomial.h
     * @author xingyifan
     * @date 2022-12-02 09:55
     *
     * @description: 分段多项式
     */

    class PiecewisePolynomial
    {
    protected:
        std::deque<Real> m_interval;   // 分段区间端点
        std::deque<Polynomial> m_poly; // 多项式列表，数量是端点数 +1，包括左右半开区间上的多项式

    public:
        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-12-02 09:58
         *
         * @description: 分段构造
         */

        explicit PiecewisePolynomial(Polynomial poly = Polynomial(0));                 // 没有区间端点，只有一个多项式
        PiecewisePolynomial(const PiecewisePolynomial &ppoly);                         // 拷贝构造
        PiecewisePolynomial(std::deque<Real> &interval);                               // 只设置区间，每段区间上是零
        PiecewisePolynomial(std::deque<Real> &interval, std::deque<Polynomial> &poly); // 设置区间和多项式

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-12-02 10:17
         *
         * @description: 重载操作符
         */

        // 多项式赋值
        void operator=(const PiecewisePolynomial &poly);

        // 多项式加减法
        PiecewisePolynomial operator+(const PiecewisePolynomial &poly) const;
        PiecewisePolynomial operator-(const PiecewisePolynomial &poly) const;
        PiecewisePolynomial operator+(const Polynomial &poly) const;
        PiecewisePolynomial operator-(const Polynomial &poly) const;

        // 多项式乘常数
        PiecewisePolynomial operator*(Real k) const;
        PiecewisePolynomial operator/(Real k) const;

        PiecewisePolynomial operator*(const PiecewisePolynomial &poly) const; // 多项式相乘
        PiecewisePolynomial operator*(const Polynomial &poly) const;          // 多项式相乘
        PiecewisePolynomial operator/(const Polynomial &poly) const;          // 多项式除法

        Real operator()(Real x) const; // 多项式求值

        // 输出函数，此函数就放在这里
        friend std::ostream &operator<<(std::ostream &o, const PiecewisePolynomial &ppoly)
        {
            for (int i = 0; i < ppoly.size() + 1; i++)
                o << ppoly[i] << std::endl;
            return o;
        }

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-12-02 10:28
         *
         * @description: 元素操作
         */

        void push_back(Real t);  // 在右端添加端点，增加一个零多项式
        void push_front(Real t); // 在左端添加端点，增加一个零多项式

        Polynomial operator[](int i) const; // 获得指定索引的多项式
        Polynomial &operator[](int i);      // 获得指定索引的多项式引用
        Polynomial poly(Real t) const;      // 获得传入点对应的多项式

        int size() const;       // 获得端点数
        int order(int i) const; // 获得指定索引的阶数

        /**
         * @file Polynomial.h
         * @author xingyifan
         * @date 2022-12-02 10:50
         *
         * @description: 多项式操作
         */

        PiecewisePolynomial intergral() const;  // 不定积分
        PiecewisePolynomial derivative() const; // 导函数

        int continuity(int i) const; // 检查端点处的连续性
        int continuity() const;      // 整体连续性

        friend PiecewisePolynomial merge(const PiecewisePolynomial &p1, const PiecewisePolynomial &p2); // 合并两个分段多项式
    };

    /**
     * @file Polynomial.h
     * @author xingyifan
     * @date 2022-12-01 22:18
     *
     * @description: 特殊的生成算法
     */

    Real factorial(int n);                                    // 阶乘
    Real combine(int n, int i);                               // 组合数
    Polynomial power(const Polynomial &poly, int n);          // 多项式整数幂
    Polynomial nodes(LinearAlgebra::Vector &c);               // 产生节点多项式
    Polynomial Lagrange(LinearAlgebra::Vector &c, int j);     // 已知节点 c0 c1 ... cn ，获得 cj 对应的初等多项式
    Polynomial Legendre(int s);                               // Gauss-Legendre 多项式
    Polynomial Bernstein(int n, int i);                       // 获得 Bernstein 基函数
    PiecewisePolynomial BasicBSpline(std::deque<Real> &knot); // 获得 B 样条基函数

    /**
     * @file Polynomial.h
     * @author xingyifan
     * @date 2022-12-09 10:34
     *
     * @description: 多项式插值算法
     */

    Polynomial interpolate(LinearAlgebra::Vector &x, LinearAlgebra::Vector &f);                                     // 多项式插值
    Real NevilleAitken(LinearAlgebra::Vector &x, LinearAlgebra::Vector &f, Real t);                                 // NevilleAitken 算法计算插值点的值
    Polynomial NewtonInterp(LinearAlgebra::Vector &x, std::vector<LinearAlgebra::Vector> &f);                       // Hermite 牛顿插值
    Table calculateTable(std::vector<Real> &x, std::vector<LinearAlgebra::Vector> &f, int width);                   // 计算差分表
    PiecewisePolynomial ppFormSpline(LinearAlgebra::Vector &x, std::vector<LinearAlgebra::Vector> &f, BCType type); // 普通 3 次样条插值
    PiecewisePolynomial cardinalBSpline2(LinearAlgebra::Vector &x, LinearAlgebra::Vector &f);                       // 2 次基数 B 样条在半整数点 1 1.5 2.5 ... 5.5 6 插值，只需要端点值
    PiecewisePolynomial cardinalBSpline3(LinearAlgebra::Vector &x, std::vector<LinearAlgebra::Vector> &f);          // 3 次基数 B 样条在整数点 1 2 3 ... 4 5 6 插值，需要两端一阶导信息
};
